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1  Technical  Summary  of  Research  Accomplished 


We  will  first  briefly  describe  the  research  results  that  have  already  been  published  in  the 
scientific  literature,  and  then  we  will  provide  some  details  of  the  remaining  publications  that 
resulted  from  the  research  effort  and  are  still  in  preparation  or  have  been  submitted. 

1.1  Summary  of  Published  Research  Accomplishments 

Paper  [1]  was  listed  in  the  report  for  Grant  F49620-99-1-0072  as  ”to  appear”  and  the 
work  presented  there  was  described  previously. 

In  paper  [2]  we  presented  an  analysis  of  the  perfectly  matched  layer  in  cylindrical  co¬ 
ordinates  discretized  with  a  staggered  second-order  accurate  finite  difference  time  domain 
method.  For  fixed  discretization  parameters,  layer  width,  and  a  quadratic  loss  function,  we 
found  that  the  numerical  reflection  produced  by  the  discrete  layer  is  accurately  predicted 
by  the  infinite  resolution  reflection  coefficient  for  amax  G  [0,  a(naJ,  where  amax  is  the  maxi¬ 
mum  value  of  the  absorption  parameter  in  the  layer.  We  also  find  that  the  finite  resolution 
reflection  coefficient  achieves  its  minimum  value  at  a  cr™ox  >  crcmax.  Numerical  experiments 
validated  the  analysis.  Figure  1  below  shows  that  the  reflection  produced  by  the  discrete 
cylindrical  PML  converges  to  zero  with  increasing  amax,  again  in  accordance  with  the  exact 
reflection  coefficient  (labelled  as  Equation  (15)  on  the  graph)=.  We  also  see  that,  for  a  fixed 
spatial  mesh  size,  A p,  the  reflection  converges  to  zero  with  increasing  layer  depth  dp,  again 
in  accordance  with  (15).  The  most  important  observation  is  the  fact  that  (15)  is  a  good 
predictor  of  the  numerical  reflection  as  a  function  of  amax  when  that  parameter  is  in  a  range 
(0,  cr^ax\.  We  expect  to  have  obtained  similar  behavior  had  we  used  dp  as  the  variable  param¬ 
eter.  A  x2  mesh  refinement  (with  all  other  parameters  fixed)  indicates  an  enlargement  of 
the  range  of  crmax  that  allows  (15)  to  be  an  accurate  predictor  of  performance.  Significantly, 
Figure  1  shows  that  Equation  22,  the  main  result  of  [2],  is  an  extremely  accurate  predictor 
of  discrete  PML  reflection. 

In  paper  [3]  we  conclusively  addressed  the  long  standing  issue  of  the  long-time  stability 
of  the  unsplit  Perfectly  Matched  Layer.  Research  groups  at  Brown  University  and  elsewhere 
had  in  the  past  noticed  linear  instabilities  in  long-time  simulations  with  th  unsplit  PML 
and  had  proposed  remedies.  Such  remedies,  however,  always  resulted  in  the  layer  not  being 
perfectly  matched.  In  [3]  we  showed  how  to  eliminate  the  undesirable  long-time  linear  growth 
of  the  electromagnetic  field  in  a  class  of  unsplit  Perfectly  Matched  Layers  (PML)  typically 
used  as  Absorbing  Boundary  Conditions  in  Computational  Electromagnetics  codes.  For 
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Figure  1:  The  maximum  of  the  actual  numerical  error  over  the  time  interval  [0, 15]  as  a 
function  of  crmax;  The  dashed  line  corresponds  to  a  x2  mesh  refinement  of  the  npml  =  10 
case.  The  straight  line  segments  represent  the  corresponding  |I?npm,|. 

the  resulting  new  PML  equations  we  gave  energy  arguments  that  showed  the  fields  in  the 
new  layer  are  bounded  by  a  time-independent  constant  hence  the  layer  is  long-time  stable. 
Numerical  experiments  confirmed  the  elimination  of  the  linear  growth,  and  the  long-time 
boundedness  of  the  fields.  For  example,  Figure  2-(a)  illustrates  the  late-time  response  of  the 
axial  Magnetic  field,  Hz,  in  the  usual  unsplit  PML;  it  is  observed  that  it  exhibits  a  linear 
growth  in  time  while  the  remaining  field  components  tend  to  zero  in  full  accordance  with  our 
analysis  in  [3].  This  simulation  was  repeated  with  7  =  0.08.  The  late-time  fields  are  shown  in 
Figure  2-(b).  Now,  all  field  components  remain  bounded  by  a  time-independent  constant  in 
full  agreement  with  our  energy  proofs  in  [3].  Further,  Figure  2-(b)  shows  that  when  7  =  0.08 
only  the  axial  magnetic  field  tends  to  a  constant  while  the  other  field  components  decay  to 
zero.  Again,  this  is  in  agreement  with  our  energy  considerations  given  in  Section  4  of  [3]. 
Significantly,  and  in  contrast  to  previous  work  on  this  topic  by  others,  our  modification  of  the 
standard  Unsplit  PML  eliminates  the  long-time  instability  while  maintaining  the  perfectly 
matched  property  of  the  resulting  layer. 

A  small  amount  of  time  was  invested  by  the  PI  during  the  academic  years  covered  by 
this  report  on  a  collaboration  with  colleagues  which  aims  to  control  fluid  interfaces  with 
electric  fields.  Papers  [7]- [8]  investigate  the  stability  of  a  thin  two-dimensional  liquid  film 
when  a  uniform  electric  field  is  applied  in  a  direction  parallel  to  the  initially  flat  bounding 
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Figure  2:  a)  The  late-time  linear  growth  of  the  axial  PML  field,  JTz(31,59,i).  The 
ifx(31,59,  t)  field  behaves  exactly  like  the  graphed  electric  field  £^(31,  59,  t).  This  case, 
7  =  0,  represents  the  standard  unsplit  PML.  b)  7  =  0.08;  the  late-time  linear  growth  of  Hz 
has  been  removed. 
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fluid  interfaces.  We  considered  the  distinct  physical  effects  of  surface  tension  and  electrically 
induced  forces  for  inviscid  and  viscous,  incompressible,  nonconducting  and  conducting,  fluids. 

1.2  Summary  of  Research  Results  to  be  Published 


In  paper  [4]  we  examine  the  short-  and  long-time  response  of  a  Cole-Cole  dielectric  half¬ 
space  subjected  to  a  delta-function  incident  pulse.  Our  purpose  is  to  contrast  the  short-  and 
long-time  impulse  response  of  the  Cole-Cole  dielectric  model  to  that  obtained  previously  by 
the  PI  for  the  Debye  model.  The  purpose  for  this  is  to  determine  whether  the  time-domain 
waveforms  obtained  in  a  Time-Domain  Reflectometry  (TDR)  experiment  could  serve  as  a 
means  to  determine  the  most  appropriate  frequency-domain  model  for  the  data  at  hand. 
Time-domain  waveforms  in  TDR  experiments  are  routinely  measured  and  analyzed.  The 
modeled  short-  and  long-time  response  of  the  experimentally  examined  dielectric  subjected 
to  any  physically  realizable  pulse  (hence  measurable  in  a  TDR  experiment)  will  be  obtained 
by  convolution  with  the  appropriate  impulse  response  function.  We  find  that  the  Cole-Cole 
impulse  response  is  infinitely  smooth  at  the  wavefront  (short-time)  in  contrast  to  the  case  of 
the  Debye  impulse  response  that  is  discontinuous  at  the  wavefront.  Also,  we  find  that  the 
location  of  the  peak  of  the  main  (late-time)  response  in  the  Cole-Cole  dielectric  occurs  at  an 
earlier  space-time  location  than  that  found  for  the  Debye  dielectric  main  response. 

In  [4]  we  formulate  the  time-domain  problem  for  the  impulse  response  of  the  Cole-Cole 
medium  and  reduce  to  a  signaling  problem  in  a  dielectric  half  space  (x  >  0)  for  the  electric 
field  which  is  given  by  the  following  Bromwich  integral 


V  ’  2'KlJC-ioo 


(1) 


where  -n  <  arg(s)  <  7 r,  ReJ-f^£  >  0,  (5  =  7  =  and  Tp  is  a  time  scale 

characterizing  the  any  other  signaling  data  which  whill  be  convolved  with  E(x,t)  above 
to  obtain  the  measured  response.  Closing  the  Bromwich  contour  with  a  semi-circle  to  the 
right  of  s  =  C  we  obtain  E(x,t)  =  0,  t  <  as  the  integrand  exhibits  no  singularities  in 
that  region.  For  the  branch  chosen,  the  square  root  in  (1)  does  not  introduce  additional 
singularities  (recall  0  <  a  <  1)  since  the  argument  of  the  roots  of  the  numerator  and 
denominator  is  arg(s)  =  (tt  +  2/c7r)/a;  k  =  0,  ±1, ...,  i.e.,  these  roots  are  outside  the  chosen 
principal  branch.  In  the  wavefront  region  (t  —*  f-)  a  large— s  expansion  of  the  bracketed 
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>  0.  Due  to0<  1  —  a  <  1,  we  determine  that 


lim  -7— E(x,t )  =  0,  n  >  0,  since  lim  sne  ^  “  =  0,  n  >  0. 
t_*e+  dtn  v  '  — 

Coo 


(2) 

(3) 

(4) 


Thus,  in  contrast  to  the  Debye  medium  model  wavefront  response,  the  response  is  infinitely 
smooth  at  the  wavefront.  For  the  special  case  a  —  1/2,  standard  Laplace  transform  tables 
allow  us  to  invert  (3)  via 


C  1{e 
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Coo  2 
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and  obtain  the  short-time  response  for  an  arbitrary  pulse  /  incident  on  the  half-space  at 
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1  0(1-7) 
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(6) 


In  fact,  by  collapsing  the  Bromwich  contour  in  (3)  onto  the  branch  cut  defined  earlier  we 
can  write  ’FQ(t)  =  £-1{e-sl  “},  t  >  0,  where 


4/Q(t)  =  -  f  e  tre  rl  D'cosf7r(1  sin  (r1  Q  sin  [7r(l  —  a)])dr.  (7) 

7T  JO 

Then,  (7)  and  the  translation  theorem  allows  us  to  write  an  integral  expression  for  E(x,  t ) 
for  any  0  <  a  <  1  as  a  convolution  of  fit)  with 

C-1{e-Aal-a}  =  A-^a(A-T^{t-— )),  t>—.  (8) 

Coo  Cqq 

For  demonstration  purposes  we  employed  Mathematica  to  evaluate  the  integral  in  (7)  for 
0  <  t  —  <  2;  we  present  in  Figures  3  and  4  the  wavefront  response  for  the  case  .4  =  1  and 

various  values  of  a;  as  a  -»  1“  the  right  hand  side  of  (7)  approaches  a  delta  function  in  time 
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Figure  3:  Wavefront  response  for  8(t)  signaling  data;  T  =  t—  -f-. 


and  the  numerical  evaluation  of  the  integral  fails.  We  verified  the  results  shown  in  Figures 
3-4  by  using  Mathematica  to  evaluate  the  following  alternative  representation  of 


oo  y.-(l~a)fc-l 

*Q(t)  =  S(“1)fc  1  T(1  +  k)T(k(a  —  1)) 


(9) 


Finally,  we  note  that  in  the  case  of  unit  step-function  signaling  data  (F(s)  =  1/s)  the  electric 
field  response  exhibits  self-similarity,  i.e,  E(x,t)  =  4?a(£),  where  £  =  A~Tr^{t  —  is  the 

similarity  variable  and  <E>a(£)  =  $  A~^'i>a(A~'^t)dt  (using  \I>Q(0)  =  0  and  the  integration 
property  of  the  Laplace  transform). 

In  the  t  >  —  region  we  evaluate  (1)  using  the  saddle  point  (steepest  descent)  method. 
We  first  rewrite  (1)  as 

E(x,t)  =  ^~  ri°°e^ds,  (10) 

where  6  =  Q(s,6 )  =  s[d  —  and  A  —  x/c ^  is  a  large  parameter.  We  obtain  the 

location  of  the  saddle  points,  s  =  s,  by  setting  8Q^6!  =  0  and  solving  for  s.  The  following 
equation  holds  for  the  saddle  points: 


Sa  +  (5  _  ap/2  Q#y/2 

sa  +  f3  7  sa  +  P  sa  +  /?7 


(11) 
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Figure  4:  Same  as  Figure  3. 

By  rationalizing  (11)  we  obtain  a  fourth  order  polynomial  in  sa  which  we  solve  for  a  repre¬ 
sentative  set  of  parameters  using  Matlab.  The  roots  are  shown  in  Figure  5.  As  9  increases 
from  1,  the  positive  root  sa  decreases  from  oo.  When  6  =  l/y^  =  Coo/c0,  sa  =  0.  For 
9  >  Cqo/co  all  real  roots  of  (11)  are  negative.  The  complex  roots  of  (11)  always  lie  in  the 
second  and  third  quadrants.  Consequently,  on  the  principal  sheet  defined  by  the  branch 
cut  introduced  earlier  by  sa,  only  the  positive  real  roots  of  (11)  survive  as  saddle  points. 
Therefore,  for  each  value  of  1  <  8  <  c^/cq  there  corresponds  one  saddle  point  oo  >  s  >  0 
on  the  real  axis.  For  9  >  CoJcq  the  saddle  point  first  coalesces  with  the  branch  point  at 
the  origin  and  then  moves  into  the  branch  cut,  hence  the  Bromwich  contour  is  no  longer 
equivalent  to  a  steepest  descent  contour.  In  the  region  1  <  9  <  Coo/co  we  apply  the  saddle 
point  method  using  X  =  x/cx  as  our  large  parameter.  A  small— s  expansion  of  Q(s,  8), 

Q{s,8)  «  s(9  -  —  +  2^3/2 s  ),  (12) 

allows  us  to  obtain  an  approximation  to  the  saddle  point,  i.e., 

s  =  (13) 


where  B  =  (1^31/^ct)  ■  We  find  that  =  a &  “  (^~  -  9)1  «  >  0  hence  the  local  steepest 
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descent  directions  at  s  are  arg(s  —  s)  =  f ,  We  obtain  the  result 

eA<3(«,0)  (1  —  a)(l  +  2a) 

£(;z,f)  «  , . ■■■■■=  -=  x  1  +  — i - V-j-  . 

\ferA <3„(J,9)  24aAB'/«(^  -  9)1+= 


In  Figure  6  we  plot  the  leading  order  term  of  (14)  as  a  function  of  a.  We  determine  that  the 
peak  of  the  long-time  response  does  not  occur  on  the  sub-characteristic  ray,  x  —  cq£,  as  in 
the  case  of  the  Debye  medium.  Also,  we  notice  that  the  leading-order  result  breaks  down, 
i.e.,  E(x,  t )  =  0,  for  6  =  Cqo/co  and  a  <  1.  This  can  be  explained  by  noticing  that  the  second 
term  in  (14),  the  correction,  diverges  at  9  =  c^/c®.  This  is  expected  as  the  derivation  of 
(14)  does  not  take  into  account  the  coalescence  of  the  saddle  point  with  the  branch  point  at 
s  —  0.  However  the  result  is  useful;  Figures  7-8  show  a  comparison  between  the  leading  order 
sadle  point  method  result  and  an  evaluation  of  (10)  for  x/c =  100  using  Mathematica. 

In  paper  [5]  we  incorporate  the  Cole-Cole  model  in  the  FD-TD  scheme.  As  a  first  step  in 
that  direction  we  developed  a  numerical  method  to  solve  the  following  fractional  differential 
equation  initial  value  problem  for  the  polarization  P 
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da  P 
dta 


+  P  —  (es  —  €oo)E;  t  >  0,  P(0)  —  0 


where 


£ P  =  1  f 1  P'W 
dt°  r(i-a)7o  (t-ty 
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Figure  6:  Leading  order  asymptotic  response  at  x/c^  =  100  as  a  function  of  a. 


is  the  Caputo  fractional  derivative  of  order  a  €  (0, 1).  We  use 


7T 


roc 

r(a)  =  /  e~zza~1dz  and  r(a)r(l  —  a)  =  — 

Jo  sm  na 

and  a  simple  change  of  variables,  z  =  £2(f  —  £),  in  the  Caputo  derivative  to  show: 

da  P  2  sin  7tq:  [°° 


roc 

Jo 


dta  7T 

where  satisfies  a  non-homogeneous  ODE  IVP 


t  >  0 


=  *  >  o,  0)  =  0;  0<£<oo. 

For  accurate  results  we  decompose  /0°°  d £  =  Jq  d£+ and  use  the  small— £  behavior  V>(£>  0  ~ 
y2a-ip^^  ancj  the  large— £  behavior  ^>(£,  t)  ~  y2Q_3P  (t),  determined  from  the  above  ODE, 
along  with  Gaussian  and  Laguerre  integration.  The  discretization  in  £  is  then  followed  by 
a  discretization  in  t  using  the  trapezoidal  rule;  the  result  is  an  0(Af2)  accurate  system 
of  JVg  +  JV[  +  1  equations  that  must  be  solved  along  with  the  FD-TD  scheme.  Typically, 
Ng  =  Ni  =  10  —  15.  We  show  tests  of  our  algorithm  for  the  problem 

da  P 

— -  +  P  =  0,  P(0+)  =  1 
dta 
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Figure  7:  Comparison  of  Equation  (14)  and  (10)  for  x/c^  =  100. 

The  exact  solution  of  the  above  equation  is 

P (t)  =  Ea(-ta),  t  >  0. 

Figures  9-10  show  very  accurate  results.  We  also  tested  the  algorithm  with  the  following 
problem 

+  P  =  H(t)  -  H(t  -  1),  P(0+)  =  0, 

where  H(t )  is  the  Heaviside  step-function.  The  picture  is  identical  to  that  in  Caputo  (Rend. 
Phys.  Acc.  Lincei,  v.  4,  pp.  89-98,  1993)  where  the  inverse  Laplace  transform  of  is 

computed,  a  =  1  corresponds  to  the  Debye  model.  Presently  [5],  a  second-order  accurate 
code  is  being  tested  that  solves  the  scaled  ID  signaling  problem,  where  u(0,t)  =  f(t), 
(x,  t)  >  0,  and  n  =  n  =  n;  =  0att  =  0 


ut  +  wt  =  vx,  vt  =  u: 


2  sin  7 ra  r00 


7T 

t2„ 


roo 

/  il;(x,^,t)d(  +  r]w(x,t)  =  nu(x,t) 

Jo 


Mx>  £>  <)  +  1 1)  =  £  t)\  0)  =  0;  0  <  C  <  oo, 


where  rj  =  ( Tp/r)a  and  k  =  77 (^  —  1). 

Finally,  a  numerical  implementation  within  a  fourth-order  accurate  FD-TD  scheme  is 
given  in  paper  [6]  for  the  the  impedance  boundary  condition  at  a  planar  interface,  separating 
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Figure  8:  Same  as  Figure  7. 


a  homogeneous  lossy  half-space  (of  conductivity  a  and  permittivity  e)  from  free-space,  which 
relates  tangential  Electric  and  Magnetic  field  components  as  follows: 

Etan  =  (J^S(t)  +  <:(t))*Htan,  (15) 

VCr 

where  rjo  is  the  impedance  of  free-space,  er  is  the  relative  permittivity  of  the  lossy  half-space, 
5(t )  is  the  delta  function,  and  ((f)  is  the  time-domain  impedance  function.  Previously  [1] 
we  have  derived  approximations  with  an  error  bound  of  Z(t)  =  +  CW)  so  the 

convolution  can  be  performed  accurately  and  quickly  since  the  approximants  are  damped 
exponentials.  Therefore,  the  remaining  obstacle  to  discretizing  (15)  to  the  order  of  accuracy 
of  the  fourth-order  scheme  is  to  correctly  extrapolate  the  tangential  magnetic  field  to  the 
boundary  node  which  is  an  electric  field  node.  For  this  purpose  we  derived  a  fifth-order 
accurate  extrapolation  formula  which  we  are  now  testing  against  a  fourth-order  extrapolation 
of  the  form  (node  0  is  the  boundary  node  where  (15)  is  to  be  implemented) 


ES+1  +  E. 


105 

=  zapprox  *  {—h1/2 

4o 


35H  4-21W 

ie"3/2  +  Ieft/2 


15 

48 


H7/2) 


where  the  r.h.s.  convolution  is  time-centered  at  f„  + 
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Figure  9:  Numerical  solution  of  +  P  =  0,  P(0+)  =  1  for  various  values  of  a. 


a=0.7 


Figure  10:  Detail  of  Figure  9  for  a  =  0.7. 
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Figure  11:  The  polarization  response  of  a  Cole-Cole  dielectric  subjected  to  a  rectangular- 
pulse  Electric  field. 

2  Conference  Papers  and  Technical  Reports 

The  following  conference  paper  was  not  reported  before  to  AFOSR,  we  include  it  here  for 
accuracy:  “Subgridding  a  Fourth-Order  FD-TD  Scheme  for  Maxwell’s  Equations,”  (with  A. 
Yefet),  Fourth  International  Workshop  Proceedings  on  Computational  Electromagnetics  in 
the  Time-Domain:  TLM,  FDTD  and  Related  Techniques,  pp.  39-45,  Nottingham,  UK,  2001. 

“Long-Time  Behavior  of  the  Unsplit  PML,  (with  E.  Becache  and  S.  D.  Gedney),  Mathemati¬ 
cal  and  Numerical  Aspects  of  Wave  Propagation  WAVES  2003,  pp.  120-124,  Springer- Verlag 
(July  2003). 

3  Presentations 

January  2002:  AFOSR  Annual  Electromagnetics  Workshop,  San  Antonio,  TX.  Title:  “A 
numerical  and  analytical  study  of  the  perfectly  matched  layer  for  Maxwell’s  equations  in 
cylindrical  coordinates.” 

February  2002:  NJIT  Applied  Mathematics  Colloquium,  Newark,  NJ.  Title:  Absorbing 
boundary  conditions  for  the  time-dependent  Maxwell  equations.” 
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June  2002:  NSF-CBMS  Regional  Conference  in  the  Mathematical  Sciences,  Numerical  Meth¬ 
ods  in  Forward  and  Inverse  Electromagnetic  Scattering,  Golden,  CO.  Title:  “A  review  of  the 
perfectly  matched  layer  ABC  and  some  new  results.” 

August  2002:  NASA  Langley  Research  Center  (ICASE),  Hampton,  VA.  Title:  “On  the 
long-time  behavior  of  unsplit  perfectly  matched  layers.” 

August  2002:  Projets  Estime  Ondes  Otto,  INRIA-Rocquencourt,  France.  Title:  “Some 
recent  results  on  unsplit  PML  for  electromagnetics.” 

January  2003:  AFOSR  Annual  Electromagnetics  Workshop,  San  Antonio,  TX.  Title:  “On 
the  long-time  behavior  of  unsplit  perfectly  matched  layers.” 

April  2003:  "Analytical  and  Computational  Methods  in  Electromagnetics,”  Spring  East¬ 
ern  Sectional  Meeting  of  the  American  Mathematical  Society,  Courant  Institute  of  Mathe¬ 
matical  Sciences,  NY.  Title:  “On  the  long-time  behavior  of  unsplit  perfectly  matched  layers.” 

June  2003:  Department  of  Mathematics,  University  of  Crete,  Heraklion,  Greece.  Title: 
“Absorbing  Layer  Boundary  Conditions  for  the  Numerical  Solution  of  the  Time-Dependent 
Maxwell  Equations  in  Open  Domains.” 

July  2003:  Sixth  International  Conference  on  Mathematical  and  Numerical  Aspects  of  Wave 
Propagation,  Jyvaskyla,  Finland.  Title:  “Long-Time  Behavior  of  the  Unsplit  PML.” 

July  2003:  Department  of  Mathematical  Sciences,  Summer  Faculty  Seminar  Series,  New 
Jersey  Institute  of  Technology,  Newark,  NJ.  Title:  “Absorbing  Boundary  Conditions  for 
Numerical  Solution  of  the  Time-Dependent  Maxwell  Equations  in  Open  Domains.” 

October  2003:  Program  in  Applied  Mathematics  Colloquia  Series,  University  of  Arizona, 
Tucson,  AZ.  Title:  Absorbing  Boundary  Conditions  for  Numerical  Solution  of  the  Time- 
Dependent  Maxwell  Equations  in  Open  Domains.” 

January  2004:  AFOSR  Annual  Electromagnetics  Workshop,  San  Antonio,  TX.  Title:  “Mod¬ 
eling  Propagation  of  Time-Domain  Pulses  in  Cole-Cole  Dielectrics.” 

April  2004:  Department  of  Mathematics  Colloquia  Series,  University  of  New  Mexico,  Al¬ 
buquerque,  NM.  Title:  “Two  Problems  in  Computational  Electromagnetics.” 
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June  2004:  Department  of  Mathematical  Sciences,  Summer  Faculty  Seminar  Series,  New 
Jersey  Institute  of  Technology,  Newark,  NJ.  Title:  “Modeling  Propagation  of  Time-Domain 
Pulses  in  Cole-Cole  Dielectrics.” 

July  2004:  EUROEM  2004  International  Conference,  Magdeburg,  Germany.  Title:  “Mod¬ 
eling  Propagation  of  Time-Domain  Pulses  in  Cole-Cole  Dielectrics.”  Invited  Presentation 
in  Special  Session  titled  ’’Time-Domain  Techniques  for  the  Transient  Analysis  of  Complex 
Problems.” 


4  Consultative,  Advisory  Functions  To  Other  Labora¬ 
tories  And  Agencies,  and  Other  Achievements 


1.  Throughout  the  period  covered  by  this  report  the  PI  continued  his  interaction  with 
Dr.  T.  M.  Roberts  of  AFRL/SNHA  on  computational  issues  that  arise  in  transient 
electromagnetic  wave  propagation.  Particularly,  I  have  been  in  communication  with 
Dr.  T.  Roberts  for  the  purpose  of  explaining/understanding  a  discrepancy  between  an 
asymptotic  (large-depth)  result,  which  predicts  decay  of  peaks  to  be  ~  (depth)-1  ^ ,  and 
an  observed  decay  of  peaks  (~  (depth)-2 /5)  obtained  in  experiments  of  electromagnetic 
pulse  propagation  in  concrete  whose  dielectric  properties  are  fitted  to  the  Lorentz 
model  of  dispersion.  For  this  collaboration  I  have  been  examining  the  wave  hierarchies 
in  Lorentz  dielectrics  exhibited  by  the  partial  differential  equation  for  the  electric  field 


e[d?(Ett  -  Ezz)  +  ccdt(Ett  -  Ezz)\  +  (Ett  -  f32Ezz)  =  0 


where  e  = 


(wg+u >l)Tl 


a 


•  and  (3  —  are  non-dimensional  parameters,  Tp  the 
incident  pulse  duration,  and  r,  u>o,  ojp,  es  are  the  Lorentz  medium  parameters  fitted  to 
experimental  dielectric  data  for  concrete  with  =  1.  For  the  particular  experimental 
setup,  e  ~  4.5  x  10~4  (small),  and  a  ~  3.3  (large)  and  /?  ~  0.6.  Using  the  method  of 
singular  perturbations  I  have  derived  the  following  partial  differential  equation  (valid 
for  |  1  and  a  «  e-1/3) 

Et  +  (3EZ  =  7  Ezz  —  5EZZZ 


where  7  and  5  depend  on  the  parameters  described  above. 
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2.  After  being  presented  with  problems  of  interest  to  Dr.  R.  E.  Peterkin,  High  Power  Mi¬ 
crowave  Division,  Air  Force  Research  Laboratory,  Kirtland  AFB  I  completed  (Novem¬ 
ber  2001)  and  forwarded  to  him  a  small  manuscript  titled  ’’Artificial  Dissipation  for 
Damping  High-Frequency  Numerical  Noise  in  FD-TD  Simulations.” 

3.  Visiting  Scientist  (no-cost  consultant),  ICASE,  NASA  Langley  Research  Center,  Hamp¬ 
ton,  VA.  July-August  2002. 

4.  Organized  a  mini-symposium  at  the  May  2004  Frontiers  in  Applied  and  Computational 
Mathematics  Conference,  entitled:  Computational  Wave  Propagation,  New  Jersey  In¬ 
stitute  of  Technology,  Newark,  NJ 

5.  Visiting  Scientist,  Projet  Ondes,  INRIA  INRIA-Rocquencourt.  Visit  funded  by  a  US- 
France  Cooperative  Research  Grant  to  Prof.  Jan  Hesthaven  (Brown  University).  July 
2004. 
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